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The equation of state and squared speed of sound c s are studied in grand canonical ensemble of 
all hadron resonances having masses < 2 GeV. The ensemble is divided into strange and non-strange 
hadron resonances. Furthermore, pionic, bosonic and femeionic sectors are considered, separately. 
It is found that c 2 calculated in the QCD matter below T c is obviously causal. There is no sign for 
C"| ■ superluminal phenomena. It is found that the lightest Goldstone bosons, the pions, represent the 

main contributors to c 2 a at low temperatures. At this temperature scale, they determine the hadronic 
thermodynamics including the equation of state, almost entirely. The comparison of the barotropic 
• dependence of the pressure calculated in the hadron resonance gas (HRG) with that of full lattice 

QCD at vanishing and finite chemical potential is excellent. Nevertheless, the comparison of c 2 a at 
vanishing chemical potential is not that good. But, when switching on the chemical potential, HRG 
c 2 and that of full lattice QCD are in good agreement, especially when c 2 a is calculated through 
the ratio s/c v , where s and c„ being entropy and specific heat, respectively. Such a discrepancy 
| in reproducing the equation of state, but not the speed of sound can be understood as the specific 

, heat reflects several types of energy susceptibilities and fluctuations. Apparently, all these collective 

■ phenomena are not entirely implemented in HRG. 

PACS numbers: 05.70.Ce, 05.45.-a, 25.75.Nq, ll.10.Wx 
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I. INTRODUCTION 



Despite the outstanding understand of the structure of matter at energy density much larger than the critical 
value which defines the hadron-quark deconfmcment phase transition, one of the yet-unsettled problems of the- 
oretical physics is the characterization of equation(s) of state (EoS) describing the behaviour of thermodynamic 
quantities at finite temperatures and densities. On one hand, the lattice QCD calculations are reliable, espe- 
cially at very high temperatures and densities. On the other hand, the exact equation of state of the hadronic 
' matter is still rather complicated. In describing the ground state properties of nuclear matter having a large 
■ number of finite nuclei, the Hartree-Fock theories using Skyrme effective interactions are shown to be quite 
successful [l|. Nevertheless, as concluded in Q, it seems that serious concerns about basic physical symmetries 
arise when using EoS derived from Skyrme interactions in framework of Hartree-Fock theories, especially at 
finite temperatures. Furthermore, it is found that the speed of sound seems to violate the causality constrains 
leading to superluminal phenomena. Regardless these difficulties in suggesting a reliable EoS for the nuclear 
ground state, the speed of sound is widely considered as an order parameter for the critical behavior in the 
strongly interacting matter. Early lattice QCD calculations [|| have shown that the speed of sound possesses 
a sharp dip over the critical region, where the deconfinemcnt phase transition is believed to take place through 
a slow crossover. Such a dip becomes weak with refining the certainty of the lattice QCD calculations 043 , 
which is based on extreme enriching the computing facilities and utilizing powerful algorithms. 

Various reasons speak for utilizing the physical resonance gas model (HRG) in predicting the hadron abun- 
dances and their thermodynamics. This model seems toprovide a good description for the thermal evolution of 
the thermodynamic quantities in the hadronic matter [9H17I] and has been successfully utilized to characterize 
the conditions deriving the chemical freeze-out at finite densities [l8l - |20j . In light of this, HRG can be used in 
calculating the speed of sound using a grand canonical partition function of an ideal gas with all experimentally 
observed states up to a certain large mass as constituents. The HRG grand canonical ensemble includes two 
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important features [121 ]; the kinetic energies and the summation over all degrees of freedom and energies of 
resonances. On other hand, it is known that the formation of resonances can only be achieved through strong 
interactions [2l| ; Resonances (fireballs) are composed of further resonances (fireballs), which in turn consist of 
resonances (fireballs) and so on. In other words, the contributions of the hadron resonances to the partition 
function are the same as that of free particles with some effective mass. At temperatures comparable to the res- 
onance half- width, the effective mass approaches the physical one [HI]. Thus, at high temperatures, the strong 
interactions are conjectured to be taken into consideration through including heavy resonances. It is found that 
hadron resonances with masses up to 2 GeV are representing suitable constituents for the partition function 
[9MT7I . Such a way, the singularity expected at the Hagedorn temperature [{J [T3| can be avoided and the strong 
interactions are assumed to be considered. Nevertheless, the validity of HRG is limited to temperatures below 
the critical one, T c . 

Recently, the problematic of characterizing a hadronic equation of state has been discussed in [22[ . In an ideal 
gas of hadron resonances with Hagedorn mass spectrum, the speed of sound is calculated at different upper 
cut-off masses in the resonance mass integration. It is found that the speed of sound initially increases similarly 
to that of an ideal pion gas, until near T c . Then, the hadron resonance effects seem to become dominant. This 
causes a vanishing speed of sound at (T c — T) 1 / 4 . 

In present work, we introduce a systematic study for the squared speed of sound based on HRG [l-[l3| ■ We 
fix cut-off resonance masses to < 2 GeV. We distinguish between the bosonic and fermionic contributions to the 
thermodynamics. Also, we distinguish between the contributions stemming from the pion gas and that from 
a gas consisting of all resonances with and without pions. Finally, we compare the results with recent lattice 
QCD calculations [ilM [23j. The causal aspects of the speed of sound are introduced in section [TTJ In a grand 
canonical ensemble, the thermal and dense evolution based on HRG calculations is derived. The HRG results 
and their comparison with full lattice QCD are given in section IIIII Section IIVI is devoted to the conclusions 
and outlook. 



II. SPEED OF SOUND IN HADRON RESONANCE GAS 



In the context of Special or General Relativity, the barotropic equation of state of a perfect fluid in the physical 
units simply reads [2J] (sf. Eq. fH) 

P = 4 e > (!) 
c A 

where the pressure p is proportional to the energy density e through the quantity oj/c 2 . This proportionality can 
be very well determined in hadronic matter and quark-gluon plasma (QGP). Furthermore, it can be determined 
over the hadron-quark phase confincment-deconfinement transition, where barotropic dependence of the pressure 
is to be deduced from the lattice QCD simulations Q. Figure Q] depicts it in a wide range of temperatures, 
1/2 < T/T c < 3. It is obvious that the dependence is almost linear referring to the nature of the phase 
transition from hadrons to quarks. This phase transition seems to be smooth i.e., simply continuous and takes 
place over a considerable range of temperature, i.e., it is slow. This kind of transitions is a very moderate than 
the second-order one, for instance. The nature of the phase diagram in lattice QCD has been discussed in (l2T |. 
In Fig. [TJ the dashed line represents the fitting in the entire T-region. Ignoring the dip around T c , the results 
can be fitted as a power law, 

p(e)= ai e a \ (2) 

where a\ = u/c 2 = 0.178 ± 0.01 and 012 = 1-12 ± 0.01, respectively. In the hadronic phase i.e., at temperatures 
< T c , the previous power law dependence seems to remain valid. Considerable changes appear in the parameters; 
ct\ = 0.096±0.003 and ct2 = 1.03±0.04. In the quark phase i.e., at temperatures > T c , the following polynomial 

p{e) = - l - + ai s a \ (3) 

describes this barotropic dependence, where a\ = 0.221 ± 0.004 and a?, = 1.072 ± 0.005. 
The squared speed of sound can be related to causality constrains [26| 

£i = J_ dp = (4) 

c 2 c 2 9(mass density) de ' 
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and therefore c s cannot exceed the speed of light c. Seeking for completeness, the propagation can be super- 
luminal if the fluid's stress tensor violates the Lorentz invariant principle (2f|. The possibility of c 2 > c 2 in 
unquantized ultradense matter has been discussed in [2t| . Furthermore, this might be the case of certain values 
for the cosmological constant, for instance, p5| . 

Apparently, expressions © and © imply that c 2 drastically changes according to the changes in the phases: 

• in hadron-quark phase: c 2 = dp/de ~ 0.199 e 012 i.e., the speed of sound varies with the thermal evolution 
of the energy density. The latter has a non-monotonic behavior, when going from hadronic to partonic 
phases and vice versa. This can be seen when comparing 

• the hadronic phase, where c 2 ~ 0.098, with 

• the partonic phase, where c 2 ~ 0.237. 

There is a drastic change in the value of c 2 , when passing through the hadron-quark phase transition. 
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Fig. 1: The pressure density, p, is drawn in dependence on the energy density e. Both quantities are given in physical 
units. Symbols are lattice QCD calculations using p4 action and temporal lattice size, N T = 8 |J. The dotted curve 
is the fitting in the quark phase, Eq. Q. The dash-dotted curve gives the fitting in the hadronic phase. The overall 
fitting is given by Eq. ([2| (double-dotted curve). The small dip at T c seems to reflect the slow phase transition known 
as crossover . 

In light of this discussion, it turns to be clear that the speed of sound seems to offer an essential tool to 
check the causality and therefore the applicability of the laws of thermodynamics, on one hand. On other 
hand, the effective equation of state of the system of interest is directly accessible through the speed of sound. 
Furthermore, it would be used as an order parameter. 

The hadron resonances treated as a free gas H3| are conjectured to add to the thermodynamic pressure 

in the hadronic phase, i.e., < T c . This statement is valid for free as well as strong interactions between the 
resonances themselves. It has been shown that the thermodynamics of strongly interacting system can also be 
approximated to an ideal gas composed of hadron resonances with masses < 2 GeV [TH, [28| . Therefore, the 
confined phase of QCD, the hadronic phase, would be modelled as a non-interacting gas of hadron resonances. 
The grand canonical partition function reads 

Z(T, V) = TV [exp- H/T ] (5) 

where H is the Hamiltonian of the system. The Hamiltonian is given by the sum of the kinetic energies of 
rclativistic Fermi and Bose particles. The main motivation of using this Hamiltonian is that it contains all 
relevant degrees of freedom of confined and strongly interacting matter. It includes implicitly the interactions 
that result in resonance formation. In addition, it has been shown that this model can submit a quite satisfactory 
description of the particle production in heavy-ion collisions. With the above assumptions the dynamics of the 
partition function can be calculated exactly and be expressed as a sum over single-particle partition functions 
Z} of all hadrons and their resonances. 



hxZ{T,n,V) = ^lnZ ? 1 (T,^) = ^±^^°°fc 2 ( ifchi{l±ex P [( All - £4 )/T]} 



(0) 



where £i(k) = (k 2 + m 2 ) 1 / 2 is the i— th particle dispersion relation, gi is spin-isospin degeneracy factor and ± 
stands for bosons and fermions, respectively. 

Before the discovery of QCD, it was speculated about a possible phase transition of a massless pion gas to a 
new phase of matter. Based on statistical models like Hagedorn [29[ and Bootstrap [3(|, the thermodynamics of 
such an ideal pion gas is studied, extensively. After the QCD, the new phase of matter is known as QGP. The 
physical picture was that at T c the additional degrees of freedom carried by QGP are to be released resulting in 
an increase in the thermodynamic quantities. The success of HRG in reproducing lattice QCD results at various 
quark flavours and masses (below T c ) changed this physical picture, drastically. Instead of releasing additional 
degrees of freedom at T > T c , it is found that the interacting system reduces its effective degrees of freedom at 
T < T c , In other words, the hadron gas has much more degrees of freedom than QGP. 

At finite temperature T and baryo-chemical potential /ij, the pressure of i-th hadron or resonance reads 



P(T, Mi) 



N a f°° 
± Y ;fVr / k 2 dk In {1 ± exp[(m 
V 2lr Jo 



(7) 



where N is the total number of hadron resonances of interest. As no phase transition is conjectured in HRG, 
summing over all hadron resonances results in the final thermodynamic pressure in the hadronic phase. Switching 
between hadron and quark chemistry is given by the correspondence between the hadronic chemical potentials 
and that of the quark constituents, for example, fii = 3nb[i q + n s fj,s, where nb(n s ) being baryon (strange) 
quantum number. The chemical potential assigned to the degenerate light quarks is fx q = (p u + /id)/2 and the 
one assigned to strange quark reads ^5 = [i q — /i s . The strangeness chemical potential ^5 is calculated as a 
function of T and /Uj under the assumption that the overall strange quantum number has to remain conserved 
in the heavy-ion collisions [l2|. Based on this assumption, /Us^/i, T) is to be calculated at each value of the 
baryo-chemical potential \i and temperature T. 

As given above, the squared speed of sound of a single boson/fermion resonance at finite T and /z is derived 
as follows. 



dp(T,(ij) f d P (T,^) \ t dT 
de{T,iXi) ~\dT) \de(T,m 



(8) 



where 
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are entropy s(T, /ij) and specific heat c v (T,Hi), respectively. The second term of Eq. (fl"U|) is complicated. Its 
physical meaning is not given, directly. But, it can be re- written as 



Cv {T,m) = ± o ^ 
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It is obvious that c v seems to reflect different types of fluctuations. The fluctuations appearing here are related 
to the energy density and its product with the chemical potential. For example, the last line gives a direct 
relation to the energy density fluctuation, \si 
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which is nothing but [(e) - Xe] T. The first line of Eq. ((TTJ) gives the averaged energy density (e) . Appendix [X] 
gives details on all these energy density dependent quantities. Therefore, it becomes obvious that the second term 
in Eq. (|TU1) indeed reflects several types of energy fluctuations. Therefore, the speed of sound calculated through 
the ratio of entropy and specific heat, Eq. (|14[) . seems to reflect essential dynamics and strong correlations 
controlling the system under investigation. Furthermore, this might interpret the results given in section IIII1 
for instance , it would explain why the QCD barotropic pressure, the equation of state, is well reproduced by 
HRG, while c 2 s is not. 



III. RESULTS AND DISCUSSIONS 
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Fig. 2: Left panel: c? s in dependence on T calculated in HRG according to Eqs. (|13|l and (|14[) . The strange hadron 
resonances are excluded. Right panel draws the same dependence but with the inclusion of strange hadron resonances. 
In both graphs, bosonic (lines) and fermionic (symbols) sectors are utilized, separately. 



The thermal evolution of the speed of sound c 2 s in natural units can approximately be given in two expressions, 
Eq. flU) and (HDJl, 

cl(T) * fg, (13) 
c 2 s(T) = (14) 

The reason of stating these two expressions is that the lattice QCD simulations [f| 0] used them. Based on 
these relations, c 2 seems to reflect not just the thermodynamics of the system. Apparently, it directly depends 
on the entropy, as well. Furthermore, c 2 is depending on the specific heat, which gives various energy density 
fluctuations, section [TTJ In Fig. EH the thermal evolution of c 2 given by Eq. (fT3")) and ([T4"|) is depicted for boson 
(dashed curves) and fcrmion (symbols) resonances, separately In the left panel, strange resonances are excluded, 
while they are included in the right panel. We find remarkable differences between bosons and fermions, between 
strange and non-strange hadron resonances and between the expressions (|13p and (|14| . It is apparent that the 
fermionic c 2 is smaller than the bosonic one. The latter decreases faster than the increase taking place in the 
earlier. The expressions (jTB")) and (fT4"|) are obviously distinguishable, especially for bosons at the intermediate 
temperatures. Including strange hadron resonances increases the fermionic c 2 , on one hand. On the other 
hand, it decreases the bosonic c 2 . In these calculations, the upper cut-off resonance mass is fixed at < 2 GeV. 
Therefore, the in/exclusion of strange resonances obviously changes the number of hadron resonances which are 
considered in the partition function, Eq. (|7|) ■ It is expected that the values of the thermodynamic quantities are 
considerably suppressed when adding new hadron resonances, i.e., additional masses, to the partition function. 
Apparently, this is not exactly the case with c 2 . An explanation would be the fact that c 2 is obtained from the 
ratio of two thermodynamic quantities, p and e or s and c v . An arbitrary change in denominator and numerator 
would not necessarily change the value of the ratio that much. 
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A. The dominant role of pions at low temperatures 




1e-05 



0.1 



0.01 



0.001 



80 



120 160 
T [MeV] 



resonances-7t 
resonances+7t 



1 20 1 60 

T [MeV] 



200 



200 



0.1 



0.01 r 



0.001 




0.0001 



80 



10 



0.1 



0.01 



80 



resonances-7i 
resonances+7: 



120 160 
T [MeV] 



200 




resonances-ji 
resonances+3i 



120 160 
T [MeV] 



200 



Fig. 3: The densities of dimensionless thermodynamic quantities p/T 4 (top left panel), e/T 4 (top right panel), s/T 3 
(bottom left panel) and c v /T 3 (bottom right) are given as functions of T for pion gas (solid curve) and resonance gas 
(dashed curves), separately. Furthermore, the pions are included and excluded from the resonance gas. The square root 
of the string tension yfa ~ 420 MeV is obtained on lattices in the continuum limit. 



Bearing in mind the effect of resonance masses on the thermodynamics, we expect that the pion mass likely 
plays the role of the effective mass at low temperatures. This has been noticed in [22j], where c 2 s of hadron 
resonance gas is found to initially increase similarly to that of an ideal pion gas. At low T, the heavy masses 
are contributing minimally to the thermodynamics. The lightest Goldstonc bosons, the pions, are dominant, 
as their massed are comparable with the temperature scale. This behavior is illustrated in Fig. [3] The four 
thermodynamic quantities mentioned so far in this section, p, e, s and c v are given in dependence on T for 
different constituents of the ideal gas described by Eq. ([7]). The units are given in y/a, which is the square root 
of the string tension yfa ~ 420 MeV obtained on lattices with certain temporal extensions and extrapolated 
to the continuum limit. Here, we distinguish between pion gas and resonance gas. The pions themselves 
can be excluded and included in the resonance gas. It is apparent that the pion gas at low T gives much 
higher thermodynamic quantities (solid curves) than that of the resonance gas (dotted curves), so that when 
including pions in the resonance gas, the values do not increase that much. Up to T ~ 120 — 140 MeV, the 
resonances without pions gives thermodynamic quantities smaller than that of pions alone. At higher T, the 
pions contributions arc no longer dominant, so that the resonance gas with or without pions comes up with 
almost the same contribution. The results of c v show that the dominance region of the pion gas is relatively 
short. 

Successfully absolving these preparations, then the thermal evolution of c 2 can be mapped out in HRG. In 
the left panel of Fig. [H we plot the same results given in the right panel of Fig. [3] In additional to the 
differentiation between bosons and fermions (symbols), both are compared with the hadrons (curves). It is 
apparent that the bosons exclusively contribute to the peak at low T. Exacter said that the pions are the main 
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Fig. 4: Left panel: c a is calculated from the ratios of p/e and s/c v in HRG and given in dependence on T for bosons 
(top symbols), fermions (bottoms symbols) and hadrons (dotted and dashed curves). A comparison between the values 
of <? a calculated in the pion (top pair of curves) and the hadron resonance gas is drawn in the right panel. When pions 
are in/excluded from the hadron resonance gas, the resulting <? a is illustrated as well (see text). 



contributors to the peak. At this temperature scale, the fermionic contributions are minimum. Increasing T 
results in a decrease in the hadronic c 2 even below the one calculated in the bosonic gas, while the fermionic c 2 
increases almost linearly slowly. 

Other features of this behavior can be revealed through the right panel of Fig. @] Here the comparison is 
made between hadronic and pionic c~. The top pair of curves shows the pionic results. It is amazing that the 
values of c 2 are the largest. The bottom pair of curves represents the results from hadron resonances excluding 
the pions. Comparing the three pions, the lightest Goldstone bosons, with the rest of hadron resonances with 
masses up to 2 GeV makes it clear that the values of c 2 are drastically reduced from ~ 0.3 for pions to ~ 0.13 for 
other resonances at T > 160 MeV. Adding the three pions to the hadron resonance results in the characteristic 
c 2 curves (middle). It starts with a peak, at low T. Increasing T causes a slow decrease in c 2 . At larger T, the 
decrease becomes fast, then becomes almost exponential. At these temperatures, the contributions of the three 
pions seem to be small, so that the results from the hadron resonances with and without pions get close to each 
others. 



B. Comparison with lattice QCD at vanishing baryo-chemical potential 

In Fig. [5j the HRG results are confronted with the full lattice QCD 0, [H|. The quark masses are given the 
physical values. We use two sets of lattice data. The first one of 2010 , where 2 + 1 staggered quark flavors and 
one-link stout improvement in the continuum limit are implemented. With 2 + 1 we mean two degenerate light 
quark flavors and one heavy strange quark. The second data set is published recently [23| . The equation of state 
of QCD for nonzero chemical potentials is determined via Taylor expansion of the pressure. Also here, 2 + 1 
quark flavors with physical masses on various lattice spacings are implemented. Although, different comparisons 
with HRG are performed for different thermodynamic quantities |23j, c 2 is not involved. 

In left panel in Fig. [5J the HRG barotropic dependence of the pressure is compared with the lattice QCD 
calculations at vanishing baryo-chemical potential. As given in Fig. [I] the pressure has an almost linear 
dependence on the energy density, especially in the hadronic phase, where the slope would be used to determine 
c 2 , directly. In the partonic phase, another relatively strong dependence seems to describe this barotropic 
relation, referring to a drastic increase in the value of c 2 . The excellent agreement between HRG and lattice QCD 
in reproducing the thermodynamic quantities does not makes it straightforward to find a proper interpretation 
for the apparent discrepancy in reproducing the speed of sound, as given in right panel of Fig. [SJ On other hand, 
it is apparent that the structure of the lattice data seems to be different than that of HRG, for instance, the 
lattice data have a minimum (dip) , while HRG does not have. A precise judgement about the data set and its 
structure are achievable through the normalized higher moments [2(|. In lattice QCD as well as in HRG, there 
is a qualitative tendency that c 2 slightly decrease with increasing T up to ~ 150 MeV. Nevertheless, HRG seems 
to overestimate the values of c 2 calculated in lattice QCD. On the other hand, there is another quantitative 
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Fig. 5: Left panel: The barotropic dependence of the pressure p as calculated in the lattice QCD 0,113] given in physical 
units (symbols) and compared with the HRG results (solid line). Right panel: in hadronic phase, the thermal evolution 
of Cg given in yfa (symbols) calculated in lattice QCD is compared with the HRG results (curves). The critical values 
are marked by the vertical lines. 



difference, this time between the HRG results themselves that are calculated by Eq. (| 13[) and Eq. (| 14[) . The 
latter is relative closer to lattice results than the earlier. 

At higher temperatures, c 2 calculated in the lattice QCD raises, while in HRG it resumes its decreasing 
tendency. In this region of temperature, where the dcconfincmcnt phase transition is assumed to set on, there is 
a clear quantitative and qualitative discrepancy The lattice results increase, so that at T ~ 185 MeV, the HRG 
results are ~ 25% smaller than that of the lattice QCD. The dip appearing in the lattice data is not present in 
HRG. The critical values characterizing the deconfinement phase transition at vanishing chemical potential are 
marked by vertical lines. 



C. Speed of sound at finite baryo-chemical potential 

The lattice QCD simulations [23| introduce a systematic estimation for c 2 in thermal and dense medium. In 
Fig. |6l we make a comparison similar to the one given in Fig. [5j but at baryo-chemical potential fib = 300 MeV. 
As discussed in section [H] the strangeness chemical potential fis has to be evaluated in dependence on fib and 
T in order to guarantee the overall strangeness conservation in the strong interactions. It is obvious that the 
speed of sound seems to be sensitive to the dense medium. Also, it is apparent that the difference between 
the two expressions (fT3|) and (fT4|) becomes large in this dense medium. It seems that c 2 resulting from the 
second expression gets suppressed when increasing baryo-chemical potential fi, while the first expression is not 
as much affected by the finite fib- This can be understood from the observation given in the left panel, where 
the barotropic relation of the pressure p, i.e., the equation of state, does not change such much when increasing 
fib from to 300 MeV. Comparing with the left panel of Fig. [5j the hadronic equation of state (solid lines) seems 
to remain almost unchanged. Despite the good agreement, the structure of the lattice data is not reproduced 
by HRG. 

It is worthwhile to notice that the HRG seems to reproduce the lattice results on c 2 at low temperatures and 
finite chemical potential. Up to T ~ 160 MeV, the agreement seems to be excellent, especially when calculating 
c 2 via the ratio s/c v . Again, at higher temperatures, a small discrepancy appears, indicating that the lattice 
results increase, while the HRG ones slightly decrease or remain constant. The discrepancy is small, because 
the critical region at this chemical potential is apparently moved to lower temperatures, ~ 160 — 165 MeV. The 
same is also valid for the lattice QCD. Furthermore, the dip observed in lattice data is positioned at lower 
temperatures, ~ 130 — 140 MeV. The critical temperature and energy density at this value of chemical potential 
are marked by the vertical dashed lines. 
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Fig. 6: The same as in Fig. but at the baryo-chemical, jib = 300 MeV. The critical values are marked by the vertical 
lines. 



IV. CONCLUSION AND OUTLOOK 



The speed of sound c 2 calculated in QCD matter below T c is obviously causal. There is no sign for superluminal 
phenomena. We conclude that the two expressions [T5l and [HI are not necessarily equivalent. Furthermore, the 
lightest Goldstone bosons, the pions, seem to represent the main contributors to the c 2 peak at low T. At this 
temperature scale, they are controlling the thermodynamics including the equation of state, almost entirely. 
The comparison of HRG with full lattice QCD at vanishing chemical potential is not satisfying. On other hand, 
when switching on the chemical potential, both results seem to have a good agreement. 

It is remarkable that the HRG can very well reproduce all thermodynamic quantities including the equation 
of state of full QCD with physical quark masses, on one hand. On the other hand, its results for c 2 seem to 
overestimate the lattice ones, especially at vanishing chemical potential. One possible interpretation would be 
the investigation of the specific heat given in Appendix [3] We show that the specific heat contains several types 
of energy density susceptibilities, fluctuations and multiplicities. It is apparent that such collective phenomena 
are not implemented in HRG, where it is guaranteed that grand canonical ensemble includes the kinetic energies 
and the summation over all degrees of freedom and energies of resonances. At high T, formation of resonances 
is only achieved through strong interactions, the effective mass is assumed to approach the physical one and 
the strong interactions are conjectured to be taken into consideration through including heavy resonances. All 
these essential aspects do not necessarily lay on collective phenomena or long range correlations. As the strong 
interaction is conjectured to be implemented, it turns to be crucial to implement other types of interactions 
and long range correlation functions. As a starting point, the effects of excluded volume and van der Waals 
repulsive interactions have to studied. In doing this, the modification from the collisionlcss, uncorrelated, ideal 
gas should be treated in framework of Uhlcnbeck-Gropper approach (3l| . Furthermore, the different energy 
density fluctuations contributing to the specific heat have to be analysed, carefully. The viscous properties in 
Hagcrdron fluid 32[ would play a crucial role in bringing HRG c 2 close to the lattice QCD one. Using different 
actions and lattice configurations, the dip appearing in the hadronic phase has to be calculated carefully in 
full lattice QCD. It is obvious that the dip has its minimum value located at a range of temperatures close to 
T c . That no dip appears in HRG would mean that the equilibrium distribution function should be a subject of 
modification [161 l33l|. Last but not least, the higher moments have to be estimated for the lattice data, in order 
to understand its structure 12011. 
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Appendix A: Energy Fluctuations in Grand Canonical Ensemble 

In a grand canonical ensemble described by the partition function ([7]), the number density 

n(T >M ) = ^Ef / fc 2 ^ , a^ , (Al) 

i=l 



and the energy density 



1 N r c 

i=l ^° 



lie t 



where is the energy of i-th state. The non-normalized second moment (known as susceptibility) |20| of this 
quantity is given by differentiating Eq. (|A2[) with respect to \i. Then 

i N roc »i- e i 1 W r oo 2"'~ e ' 

*^)=^£f / ^-^^T^Ef/ *^ 7 6 (A3) 

= ( e ) T 'y^ / k 2 dk ei — 

i=i J ° (l±e 



e T 

a, (A4) 



where the second term in Eq. (|A4[) is exactly the last line in Eq. (fTTj) multiplied by T. 

Furthermore, the three integrals appearing in the second line of Eq. (jlljl have a common fraction, so that 
this line can be approximately summarized as follows. 



k 2 dk O 5 , (A5) 

(l±eV 



where E [— £i,s? — SiHi]- A similar integral seems to appear in the product of the susceptibility of O and the 
number density [20[, 



Xo no 



^ 27rT 2 7 1 / 1±e (^- £l )/n 2 [27rT7 



2ttT 7 lie( £ «-M0/r 



(A6) 



The question arises now is: "What is the product of the susceptibility of O and the number density?" This can 
be found through the example 

X n = [(AO - (TV 2 )] (AO (A7) 

where the multiplication of two averages is obviously commutative. Generally it results in summation of the 
product permutation of the two sequences times the inverse product of their length. Therefore, 

X n = (N) 2 -(N 2 )(N). (A8) 



It is apparent that each of these three terms, Eq. (|A5|) . gives a certain type of energy fluctuations or multiplic- 
ities, while the fourth term given Eqs. (|A4[) and (fT2j) is obviously directly related to the energy susceptibility. 
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